Probabilistic Design Tool for Optimizing a Technical System

ABSTRACT

A nonlinear technical product or process described by stochastic system output target values dependent on stochastic system input parameter values, thereby stating discrete technical system dependencies, is optimized by using a Response Surface Methods based on discrete technical system dependencies to generate at least one continuous auxiliary function for the real dependencies of the target values on the input parameter values. Next, an auxiliary function is used to generate at least one optimizing parameter for optimization by an objective function, thereby generating an additional discrete technical system dependence. The technical system is adaptively optimized by repeating the above, using the additional discrete technical system dependence, until the difference of successive optimized optimizing parameters is below a threshold. The final additional discrete technical system dependence is an optimal technical system operating point.

TECHNICAL FIELD

The present application relates to a method according to introductionpart of the main claim and to a corresponding application.

STATE OF THE ART

It is already well known that tolerances of important influence factorsmust be taken into account for planning a technical system or atechnical product. Conventionally high security margins are provided totake no risk respectively to take risks as low as possible for theplanning of a system. This may lead to high fabrication and/or operationcosts.

For some special applications conventional software packages do exist:

COMREL

COMREL is based is based on FORM/SORM and exists in two variants(FORM/SORM are first and second order reliability methods). COMREL isfor reliability analysis of components. COMREL consists of two parts:COMREL-TI for time invariant and COMREL-TV for time variant reliabilityanalysis. Base for both program parts is the method of first order(FORM) or second order (SORM). COMREL-TI can be supplied separately.COMREL-TV bases on COMREL-TI. COMREL comprises two alternative,efficient and robust algorithms to find the so called beta-point (pointfor the local highest constraint or failure probability). Thebetter-point is the base for the FORM/SORM method for a probabilityintegration. Other options for probability integration are mean valuefirst order (MVFO), Monte Carlo simulation, adaptive simulation, sphericsimulation and several importance sampling schemata. 44 differentprobability distributions (univariate stochastic models) are useable.Arbitrary dependency structures can be generated with the aid of theRosenblatt-transformation, the equivalent correlation coefficientsaccording to Nataf, Der Kiuerghian or of the Hermite-models. Next to thereliability index importance values for all relevant input values arecalculated: Global influence of the basis variables to the reliabilityindex, sensitivities and elasticities, for the distribution parameters,the mean values and the standard deviations of the basis variables,sensitivities and elasticities for deterministic parameters in theconstraint or failure function. Out of the sensitivity analysis partialsecurity factors are deviated. Parameter studies can be performed forarbitrary values, e.g. for a distribution parameter, a correlationcoefficient or a deterministic parameter. Basing on a parameter studycharts of the reliability index, of the failure or of the survivalprobability, of the influences of basic variables or of deterministicparameters and of the expectancy value of a cost function can begenerated. All results are available as a structured text file and as afile for the generation of plots. Charts generated in COMREL andformatted with the extensive plot options can be exported with the usualWindows equipment (Clipboard, Metafile, Bitmap). If it is necessary, adetailed outprint of provisional results for a failure search can begenerated.

NESSUS

NESSUS (Numerical Evaluation of Stochastic Structures Under Stress) isan integrated finite element program with probabilistic load.Probabilistic load. Probabilistic sensitivities relating to μ and σ:FORM/SORM/FPI (fast probability integration), . . . connection to ANSYS,ABAQUS, DYNA3D. NESSUS is a modular computer software system forperforming probabilistic analysis of structural/mechanical componentsand systems. NESSUS combines state of the art probabilistic algorithmswith general-purpose numerical analysis methods to compute theprobabilistic response and reliability of engineered systems.Uncertainty in loading, material properties, geometry, boundaryconditions and initial conditions can be simulated. Many deterministicmodeling tools can be used such as finite element, boundary element,hydrocodes, and user-defined Fortran subroutines.

DARWIN

DARWIN (Design Assessment of Reliability With Inspection) This softwareintegrates finite element stress analysis results, fracture-mechanicsbased life assessment for low-cycle fatigue, material anomaly data,probability of anomaly detection and inspection schedules to determinethe probability of fracture of a rotor disc as a function of appliedoperating cycles. The program also indicates the regions of the diskmost likely to fail, and the risk reduction associated with single andmultiple inspections. This software will be enhanced to handle anomaliesin cast/wrought and powder nickel disks and manufacturing andmaintenance-induced surface defects in all disk materials in the nearfuture.

The programs NESSUS, DARWIN and COMREL have a certain distribution inindustry. All those programs merely concern mechanical reliabilityanalyze. Finite element packages are integrated, in which stochastic isdirectly integrated. Thus the stochastic distribution of the load candirectly be converted into the distributions of the displacements and acomponent part can be converted into risk zones by this. For elected,external finite element programs there exist interfaces at NESSUS andCOMREL. DARWIN and COMREL moreover merely offer an instationaryanalysis. The process variables can be stochastic processes on a limitedscale. Stochastic optimization, is not integrated within NESSUS, DARWINand COMREL.

Nonlinear optimization algorithms cope with the problem:

$\begin{matrix}{{\min\limits_{\overset{->}{x} \in ^{n}}{f\left( \overset{->}{x} \right)}},{{g\left( \overset{->}{x} \right)}\underset{\_}{>}0},} & (1)\end{matrix}$

with g({right arrow over (x)}) is a constraint, especially a failure, ofan arbitrary value, the constraint or failure being caused bydeterministic input parameters.

When {right arrow over (x)} are no longer deterministic variables butstochastic random variables (e.g. normal distributed random variable{right arrow over (x)}ε N({right arrow over (μ)}, Σ)) the deterministicoptimization problem (1) passes into the following probabilisticoptimization problem:

$\begin{matrix}{{\min\limits_{\overset{->}{x} \in {N{({\overset{->}{\mu},\Sigma})}}}{E\left( {f\left( \overset{->}{x} \right)} \right)}},{{P\left( {{g\left( \overset{->}{x} \right)}\underset{\_}{<}0} \right)}\underset{\_}{<}{tol}}} & (2)\end{matrix}$

I.e. the expectation value respectively mean value of the target size,E(f({right arrow over (x)})), is minimized and the constraints may beviolated up to a prescribed tol. The mean values {right arrow over (μ)}of the input parameters are the design parameters. FIG. 1 illustratesthe optimization problems (1) and (2). The deterministic optimizationproblem passes into the probabilistic optimization problem. When designvariables respectively input parameters are afflicted with uncertaintiesthe optional dimensioning of a technical system y=f({right arrow over(x)}) leads to a new and different operating point. In other wordsminimizing the function f({right arrow over (x)}) leads to results thatare different from those obtained by minimizing the mean value of thefunction f({right arrow over (x)}). This effect may be observed in FIG.2. FIG. 2 shows an optimal dimensioning of the system y=f({right arrowover (x)}), whereby an deterministic optimization is shown left, arepresentation of the constraint or the failure in case of adeterministic optimization is shown in the middle and a probabilisticoptimization of the system is shown on the right.

A popular method for the computation of the response of a stochasticsystem is the Monte-Carlo method. The computation of the mean value andthe variance of the system y=f({right arrow over (x)}) is presented inthe following Table:

Monte Carlo Method:

START: Determine a set {right arrow over (x)}₁, . . . , {right arrowover (x)}_(m), which represents the distribution of the input parameter.

$\begin{matrix}{{{{{ITERATE}\text{:}\mspace{20mu} j} = 1},2,3,...\mspace{14mu},m}\text{}{y_{j} = {f\left( {\overset{->}{x}}_{j} \right)}}{{END}\mspace{20mu} j}\text{}{{{EVALUATE}\text{:}\mspace{20mu} {E(Y)}} = {\frac{1}{m}{\sum\limits_{1\underset{\_}{<}j\underset{\_}{<}m}y_{j}}}}{{{EVALUATE}\text{:}\mspace{20mu} {V(Y)}} = {\frac{1}{m - 1}{\sum\limits_{1\underset{\_}{<}j\underset{\_}{<}m}\left( {y_{j} - {E(Y)}} \right)^{2}}}}} & \;\end{matrix}$

In order to assure a correct computation of the characteristic sizesE(Y) and V(Y), the size m of the ensemble must be very large. Hence,embedding the Monte-Carlo method into a framework of optimization isdifficult in practical cases: To handle computational fluid dynamics orlarge Finite Element problems in reasonable time neither a supercomputer nor a large cluster of workstations would suffice.

DISCLOSURE OF THE INVENTION

It is an object of the present invention to efficiently reduce costs fordesigning nonlinear technical systems, especially like technicalproducts or technical processes. Especially optimized operating pointsof the technical systems should be easily and efficiently found in ashort period of time. Especially computational fluid dynamics or largeFinite Element problems should be handled in reasonable time.Accordingly the optimization of the system should be “time efficient”,that is the necessary period of time for achieving an optimized resultshould be short in comparison with conventional methods, for examplewith the Monte-Carlo method.

The object of the present invention is solved by a method according tothe main claim. Advantageous embodiments are claimed by the subclaims.

The present method was developed in order to optimize nonlineartechnical systems which are afflicted with uncertainties. Inputparameters or model parameters of general technical systems mayfluctuate, i.e. may have an underlying stochastic distribution. Theseuncertainties of the input parameters are carried over to the targetvalues, which also have underlying stochastic distributions. Acontinuous auxiliary function is calculated on the base of thesestochastic dependencies. Afterwards using the auxiliary function ({tildeover (f)}_(i)({right arrow over (x)}); {tilde over (g)}_(i)({right arrowover (x)})) at least one optimizing parameter (E({tilde over(f)}_(i)({right arrow over (x)})); Var({tilde over (f)}_(i)({right arrowover (x)})); P({tilde over (g)}_(i)({right arrow over (x)}))) isgenerated. An objective function is used for optimizing the optimizingparameter thereby generating a discrete technical system dependence({tilde over (f)}_(i)({right arrow over (x)}_(i)); {tilde over(g)}({right arrow over (x)}_(i))). This dependence corresponds to aninterpolation point. The newly generated interpolation point is used formaking the stochastic dependencies more accurate by adding theinterpolation point to the stochastic dependencies of the technicalsystem. Again a continuous auxiliary function is calculated byinterpolation to repeat the two step cycle. The cycle can be repeateduntil the difference of successive optimized optimizing parameters(|E({tilde over (f)}_(i)({right arrow over (x)}_(i))−E(f_(i+1)({rightarrow over (x)}_(i+1))|; |Var({tilde over (f)}_(i)({right arrow over(x)})−Var({tilde over (f)}_(i+1)({right arrow over (x)}_(i+1))|;|P({tilde over (g)}_(i)({right arrow over (x)}_(i)))−P({tilde over(g)}_(i+1)({right arrow over (x)}_(i+1)))|) is as low as desired. Athereto belonging last additional discrete technical system dependence({tilde over (f)}_(p)({right arrow over (x)}_(p)); {tilde over(g)}_(p)({right arrow over (x)}_(p))) is useable as an optimal technicalsystem operating point. Thus the technical system described by certainphysical values (physical values of the nonlinear technical system canbe length, area, volume, angle, time, frequency, velocity, acceleration,mass, density, force, moment, work, energy, power, pressure, tension,viscosity, and all further physical kinds of quantities) is optimizedstochastically. For example {tilde over (f)}_(p)({right arrow over(x)}_(p)) is the transmitting power of a transmitter depending on thearea and/or the align angle of an antenna. Knowing the method accordingto the present invention the one skilled in the art is able tostochastically optimize arbitrary technical systems basing on technicaland/or physical parameters without being inventive.

Alternative or cumulative objective functions are claimed. On thecontrary to the state of the art the present method is a stochasticoptimizer allowing a common interface. In this case the present methoddoes optimization without modeling.

The focal point clearly lies on stochastic optimization. An additionaldiscrete technical system dependence ({tilde over (f)}_(i)({right arrowover (x)}_(i)); {tilde over (g)}_(i)({right arrow over (x)}_(i))) can bea base for an additional interpolation point being used for calculatinga continuous auxiliary function with a higher accuracy than theprecedent continuous auxiliary function. Repetition is performed byexecuting a first two step cycle (i=1) followed by a second two stepcycle (i=2) followed by third two step cycle (i=3) and so on up to alast two step cycle with i=p. Accordingly i=1, 2, 3 . . . p or in otherterms {iεN/1, 2 . . . p}.

According to the present method an improved approach is proposed. Basingon the present method optimizing the mean value and optimizing thevariance of the target parameter are two efficient alternatives forobjective functions. The objective function can be provided byoptimizing the mean value (E(f({right arrow over (x)})) determined bythe formula

E(y)=

(f({right arrow over (x)})ρ({right arrow over (x)})d{right arrow over(x)})  (3)

and/or by optimizing the variance (Var(f({right arrow over (x)})) of thetarget value (y=f({right arrow over (x)})) determined by the formula

Var(y)=

(f({right arrow over (x)})−E(y))²ρ({right arrow over (x)})d{right arrowover (x)},  (4)

where ρ({right arrow over (x)}) is a probabilistic density function ofthe input parameter distribution. The precedent integrals for the systemy=f({right arrow over (x)}) are calculated numerically. The efficientcalculation of the mean value (expectation value) and of the variance ofa system is discussed in the DE 10308314.6 “Statistische Analyse einestechnischen Ausgangsparameters unter Berücksichtigung der Sensitivität”the content of which is hereby totally introduced into the disclosure ofthe present description.

According to an advantageous embodiment optimizing is performed byminimizing or maximizing the mean value (E(f({right arrow over (x)}))and/or minimizing the amount of the variance (|Var(f({right arrow over(x)})|) of the target value (y=f({right arrow over (x)})).

According to another advantageous embodiment the objective function isalternatively or cumulatively provided by optimizing the input parameter({right arrow over (x)}) tolerances (σ_(i)=1, . . . , n), with the inputparameter tolerances (σ_(i)) are especially maximized

$\begin{matrix}{\left( {\max\limits_{\sigma_{i},{i = 1},\ldots \mspace{14mu},n}{\sum\limits_{i = 1}^{n}\sigma_{i}}} \right).} & (5)\end{matrix}$

The calculation of the maximally allowed tolerances with the presentmethod can provide efficient reductions of costs.

According to another advantageous embodiment optimizing is performed byadditionally keeping a constraint like a failure probability

(P(g({right arrow over (x)})≦0)=∫_(g({right arrow over (x)})≦0)ρ({rightarrow over (x)})d{right arrow over (x)}))  (6)

of another value under or equal to a prescribed probability tolerance(tol), with ρ({right arrow over (x)}) is a probabilistic densityfunction of the input parameter distribution. Thus general nonlineardeterministic constraints may be added to the optimization problem.Accordingly the present method is capable to treat stochasticconstraints. Hence, the present method can keep e.g. a constraintespecially the failure probability in a mechanical system within aprescribed limit of tolerance. Keeping a constraint or failureprobability under or equal to a prescribed probability tolerance (tol)means keeping a probability of differences to constraints (P(g({rightarrow over (x)})≦0)) under or equal to a prescribed probabilitytolerance (tol). Constraints can be nonlinear deterministic and/orstochastic.

According to another advantageous embodiment the objective function isalternatively or cumulatively provided by optimizing, especiallyminimizing, a constraint or failure probability (P(g({right arrow over(x)})≦0)=∫_(g({right arrow over (x)})≦0)ρ({right arrow over (x)})d{rightarrow over (x)}) of another value, with ρ({right arrow over (x)}) is aprobabilistic density function of the input parameter distribution

According to the advantageous embodiment of mixed stochastic inputparameters ({right arrow over (x)}) and deterministic input parameters({right arrow over (x)}_(D)) in the formulas ({right arrow over (x)}) issubstituted to ({right arrow over (x)},{right arrow over (x)}_(D))and/or a constraint like the failure of another value caused by thedeterministic input parameters (h({right arrow over (x)}_(D)) is limitedto ≦0.

According to another advantageous embodiment optimizing is performed byusing sensitivities

$\begin{matrix}{\frac{\partial f}{\partial x_{i}},{i = 1},...\mspace{14mu},n} & (7)\end{matrix}$

of the target value (y=f({right arrow over (x)})) with respect to theinput parameters ({right arrow over (x)}). As for the deterministicoptimization (1), the sensitivities with respect to the input variablesare required. The efficient calculation of the mean value (expectationvalue) and of the variance of a system is discussed in the DE 10308314.6“Statistische Analyse eines technischen Ausgangsparameters unterBerücksichtigung der Sensitivität” the content of which is again herebytotally introduced into the disclosure of the present description.

The present probabilistic optimizer method is able to treat withstochastic design variables with normal distributions and betadistributions for the design variables. Both distributions may behandled at the same time and the normal distributed variables may alsobe dependent. The normal distributed design variables have the density:

$\begin{matrix}{{\rho \left( \overset{->}{x} \right)} = {\frac{\sqrt{{Det}\; \Sigma}}{\left( \sqrt{2\pi} \right)^{n}}{\exp \left( {{- \frac{1}{2}}\left( {\overset{->}{x} - \overset{->}{\mu}} \right)^{T}{\sum^{- 1}\left( {\overset{->}{x} - \overset{->}{\mu}} \right)}} \right)}}} & (8)\end{matrix}$

The beta distributed design variables have the density:

$\begin{matrix}{{\rho \left( \overset{->}{x} \right)} = {\prod_{i = 1}^{n}\; {\frac{\Gamma \left( {\alpha_{i} + \beta_{i}} \right)}{{\Gamma \left( \alpha_{i} \right)}{\Gamma\beta}_{i}}{x_{i}^{{\alpha \; i} - 1}\left( {1 - x_{i}} \right)}^{\beta_{i} - 1}}}} & (9)\end{matrix}$

The beta distribution has the advantage that also asymmetricdistributions and as a special case even the uniform distribution can berepresented. If the input distributions are given in terms of discretepoints, the parameters of the normal distribution may be identified bythe Gauss-Newton algorithm.

According to another advantageous embodiment firstly a stochasticevaluating of the technical system is performed on the base of anonlinear technical system function and of the density functions of thesystem input parameters

${\rho \left( \overset{->}{x} \right)} = {\frac{\sqrt{{Det}\; \Sigma}}{\left( \sqrt{2\pi} \right)^{n}}{\exp \left( {{- \frac{1}{2}}\left( {\overset{->}{x} - \overset{->}{\mu}} \right)^{T}{\Sigma^{- 1}\left( {\overset{->}{x} - \overset{->}{\mu}} \right)}} \right)}}$

and/or

${\rho \left( \overset{->}{x} \right)} = {\prod_{i = 1}^{n}{\frac{\Gamma \left( {\alpha_{i} + \beta_{i}} \right)}{{\Gamma \left( \alpha_{i} \right)}{\Gamma\beta}_{i}}{x_{i}^{{\alpha \; i} - 1}\left( {1 - x_{i}} \right)}^{\beta_{i} - 1}}}$

by calculating the mean value E(y)=

(f({right arrow over (x)})ρ({right arrow over (x)}), the varianceVar(y)=

(f({right arrow over (x)})−E(y))²ρ({right arrow over (x)})d{right arrowover (x)}, the density function and the cumulative density function ofthe technical system response. Different to deterministic analysis, herethe sensitivities are required not only for the stochastic optimization,but also for the stochastic evaluation.

According to another advantageous embodiment in case of discretedistributed input parameters ({right arrow over (x)}_(j) with j=1, 2, 3,. . . , m) and corresponding discrete target values (y_(j)=f({rightarrow over (x)}_(j))), the following steps can be performed bygenerating of an nonlinear auxiliary model for the technical product orthe technical process, especially by an polynomial approximation to thediscrete data; further by time efficient optimizing the technical systemusing one or more of the objective functions comprising a statisticrepresentation, thereby generating an operating point of the nonlineartechnical system. Especially “high dimensional” input parameters ({rightarrow over (x)}_(j)), with j=1, 2, 3, . . . , m and for example withm≧20 (for “j” and “m” see the above described Monte-Carlo method) shouldbe handled relatively rapidly for example in comparison with theMonte-Carlo method. Other possibilities are for example m≧30, m≧50,m≧100, m≧500, m≧1000 . . . .

According to another advantageous embodiment the input parameters({right arrow over (x)}) satisfy common stochastical differentialequations, whose density development is described by the Fokker-Plankequation.

According to further advantageous embodiments the Response SurfaceMethods are the Hermite, the Laguerre, the Legendre and/or the Jacobiapproximation. Other conventional methods are also possible.

The present invention is described in connection with embodiments whichare shown in the Figs. They show:

FIG. 1 a comparison of the conventional deterministic system with theprobabilistic system;

FIG. 2 a comparison of the conventional deterministic optimizationproblem with the claimed probabilistic optimization problem;

FIG. 3 a comparison of the conventional Monte-Carlo evaluation with theapproximations according to the present method;

FIG. 4 a diagram of a pareto set of two objectives;

FIG. 5 a second order polynomial approximation to discrete stochasticdata;

The difference of the present method approximates and the Monte-Carloevaluations may be observed in FIG. 3.

Probabilistic Design Goals

The deterministic optimization problem (1) determines the optimaloperation point. When changing the probabilistic design parameters someSix Sigma relevant goals respectively probabilistic design goals can beachieved which could be treated by the present probabilistic optimizermethod. According to the used objective function(s) the followingembodiments for optimization goals can be achieved:

1. Stochastic evaluation of nonlinear systems: Given a generalnon-linear system and the density functions of the system inputparameters, the present method is able to compute the stochasticresponse of the system without any Monte-Carlo evaluation. To beconcrete, the mean value, the variance, the density function and thecumulative density function of the system response may be computed.Thus, parametric studies of the system may be performed.2. Probabilistic optimization: The mean value of the system y=f({rightarrow over (x)}) is minimized. Additionally, another value could belimited to a given probability (constraint or failure probability).

$\begin{matrix}{{\min\limits_{\overset{->}{x} \in {N{({\overset{->}{\mu},\Sigma})}}}{E\left( {f\left( \overset{->}{x} \right)} \right)}},{{P\left( {{g\left( \overset{->}{x} \right)}\underset{\_}{<}0} \right)}\underset{\_}{<}{tol}}} & (10)\end{matrix}$

3. Robust design: The variance of the system y=f({right arrow over (x)})is minimized. That is, the system is shifted into states, which are notsensitive with respect to perturbations of the input parameters.Additionally, another value could be limited to a given probability(constraint or failure probability).

$\begin{matrix}{{\min\limits_{\overset{->}{x} \in {N{({\overset{->}{\mu},\Sigma})}}}{{Var}\left( {f\left( \overset{->}{x} \right)} \right)}},{{P\left( {{g\left( \overset{->}{x} \right)}\underset{\_}{<}0} \right)}\underset{\_}{<}{tol}}} & (11)\end{matrix}$

4. Robust design optimization: Any combination of the preceding casesmay be optimized:

$\begin{matrix}{{{\min\limits_{\overset{->}{x} \in \; {N{({\overset{->}{\mu},\Sigma})}}}{\alpha \; {E\left( {f\left( \overset{->}{x} \right)} \right)}}} + {\beta \; {{Var}\left( {f\left( \overset{->}{x} \right)} \right)}}},{{P\left( {{g\left( \overset{->}{x} \right)}\underset{\_}{<}0} \right)}\underset{\_}{<}{tol}}} & (12)\end{matrix}$

The variables “α” and “β” are weighting factors for weighting the meanvalue and the variance. Looking for a robust operation point and lookingfor a probabilistic optimal operating point may be competing targets.Therefore, all combinations of the weighted sum (12) may be reasonable.In a further embodiment of the present method, the Pareto set of the twoobjectives is computed, see FIG. 4. FIG. 4 shows the computing thepareto set of probabilistic optimality (maxE(f({right arrow over (x)})))and robustness (max−Var(f({right arrow over (x)}))).

5. Minimization of constraint or failure probability: In many cases, itmakes sense to minimize the failure probability directly (instead oflimiting failure probability by a given value).

$\begin{matrix}{\min\limits_{\overset{->}{x} \in \; {N{({\overset{->}{\mu},\Sigma})}}}{P\left( {{g\left( \overset{->}{x} \right)}\underset{\_}{<}0} \right)}} & (13)\end{matrix}$

Constraints can be also optimized by maximizing a constraintprobability.

6. Maximization of input tolerances: Questions of cost lead to thefollowing problem: How inaccurate may a system or product be producedwhile keeping its constraint or failure probability within a giventolerance? Let {right arrow over (x)} be independent random variables,e.g. {right arrow over (x)}εN({right arrow over (μ)}diag(σ₁, . . .σ_(n))).

The question is how large the variances σ_(i)=1, . . . , n may be chosenwhile keeping the constraint or failure probability within a giventolerance.

$\begin{matrix}{{\max\limits_{\sigma_{i},{i = 1},\ldots,n}{\sum\limits_{i = 1}^{n}\sigma_{i}}},{{P\left( {{g\left( \overset{\rightarrow}{x} \right)} \leq 0} \right)} \leq {tol}}} & (14)\end{matrix}$

7. Mixed deterministic and probabilistic Design variables: Modelingtechnical systems often deterministic and probabilistic design variablesarise at the same time. The optimization problem (11) then becomes:

$\begin{matrix}\begin{matrix}\; & {{\alpha \; {E\left( {f\left( {\overset{\rightarrow}{x},{\overset{\rightarrow}{x}}_{D}} \right)} \right)}} + {{\beta Var}\left( {f\left( {\overset{\rightarrow}{x},{\overset{\rightarrow}{x}}_{D}} \right)} \right)}} \\\min\limits_{{({x_{D},\overset{\rightarrow}{x}})} \in {({{N{({\overset{\rightarrow}{\mu},\Sigma})}},^{''}})}} & {{P\left( {{g\left( {\overset{\rightarrow}{x},{\overset{\rightarrow}{x}}_{D}} \right)} \leq 0} \right)} \leq {tol}} \\\; & {{h\left( {\overset{\rightarrow}{x}}_{D} \right)} \leq 0}\end{matrix} & (15)\end{matrix}$

The present method is able to treat the above optimization problem (15).The variables “α” and “β” are weighting factors for weighting the meanvalue and the variance. The density function and the accumulated densityfunction of system output y=f({right arrow over (x)}) are calculatednumerically in every case. The stochastic sensitivities, that is, thederivatives of the output moments with respect to the input moments area byproduct of optimization. They are available in every state of thesystem. Point (1)-(7) suggest, that the system input must be normaldistributed ({right arrow over (x)})δ(N({right arrow over (μ)}, Σ)). Thepresent method is able to treat also with mixed normal and betadistributions.

Highlights of the present method are especially:

-   -   the multi critical optimization of robustness and stochastical        optimality.    -   and the maximization of tolerance of the process parameters        respectively of the input parameters at predetermined limiting        of the reliability probability.

Basing on the present method also instationary analysis should beperformed. Therewith the process variable respectively the inputvariable can satisfy common stochastical differential equations, whosedensity development is described by the Fokker-Plank equation. Aninstationary optimization, e.g. the optimization of the period of life,is not known by the state of the art.

Design of Experiments (DOE)

It may happen, that no physical model is available for a complicatedprocess. In this case, the present method is able to construct aauxiliary model from discrete date of the system. With this auxiliarymodel, all the analysis of the present method, given in the last sectionmay be performed. Of course the validity of such a model is only givenin a small range. To demonstrate this, a comparison of the nonlinearmodel with the auxiliary model is given.

Consider a very simple nonlinear model is given by

f(x,y)=(exp(−3*x)+2*arc tan(x)+exp(8−c))*(y*y+1)  (16)

X and Y are normal distributed random variables with

$\begin{matrix}{{\mu_{1} = {{E(X)} = 7.0}},{\mu_{2} = {{E(Y)} = 0.5}},{\Sigma = \begin{pmatrix}0.5 & 0.01 \\0.01 & 0.1\end{pmatrix}}} & (17)\end{matrix}$

A stochastic analysis by RODEO gives the corresponding mean value of thenonlinear system:

E(f(x,y))=4.47  (18)

Now we want to minimize the system f(x,y) in the stochastical sense.

In a first step a deterministic optimization is performed:

$\begin{matrix}{\min\limits_{x,y}{f\left( {x,y} \right)}} & (19)\end{matrix}$

leads to the values

x=0.14,y=0.0,f(x,y)=0.93,E(f(x,y))=2.37  (20)

Using stochastic optimization by the present method, see last sections,

$\begin{matrix}{\min\limits_{\mu_{1},\mu_{2}}{E\left( {f\left( {x,y} \right)} \right)}} & (21)\end{matrix}$

leads to the values

μ₁=0.6,μ=0.0,f(x,y)=1.25,E(f(x,y))=1.42  (22)

First we could state that stochastic optimization results in a higherdeterministic value (f(x,y)=1.25) but in a much smaller stochastic valueE(f(x,y))=1.42).

In the next step, we assume that we have no longer a nonlinear model butdiscrete normal distributed values (17) and additionally the discretecorresponding system response. With a random generator normaldistributed values are generated in the range of (x,y)ε[5:95]×[−0.4;1,4]. The present method is able to fit an auxiliary model to this data.FIG. 5 shows a second order polynomial approximation of the data. FIG. 5presents a nonlinear model according to (16) versus a fitted auxiliarymodel from discrete data. Upper diagram shows that the quadricapproximation becomes bad for values x≦5. Middle diagram shows a goodapproximation in the range (x,y) are elements of [5:9.5]x[−0.4:1.4]. Thelower diagram shows a zoom into the range (x,y) are elements of[5:9.5]x[−0.4:1.4].

Also with this model, the present method is able to improve theoperating point.

Stochastic minimization of the mean value leads to

μ₁=5.7,μ_(2=0.0) ,E(f(x,y))=3.2  (23)

which is an improvement to (18).

Applications of the Present Method

To sum up Design For Six Sigma (DFSS) or probabilistic design is a taskof ongoing interest when manufacturing products or controllingprocesses. These methodologies try to analyze in which way uncertaintiesof the design parameter influence the system response and try to improvethe system or product. The present probabilistic optimizer method isdesigned to support the goals of Six Sigma, see section “probabilisticdesign goals”. There are two main applications of the present method:

-   -   Many technical processes in the field of aerodynamics,        electromagnetics and structural mechanics can be simulated by        software packages, which may be distributed commercially.        Depending on the discretization and the method itself much        computation time is used to solve these problems. The present        method is designed to treat with these problems by avoiding        expensive Monte-Carlo evaluations.    -   For many processes no physical model is available. Only        measurements of the design parameter and the process response        define the process. The present method is able to construct an        auxiliary model by using polynomial approximation (DOE). With        this pseudo-model all the analysis given in section        “probabilistic design goals” may be performed.

Generally the present method was developed to optimize systems orproducts whose influence parameters fluctuate. Optimizing can mean, thatthe system or product is provided as robust as possible or are providedas optimal as possible in a probabilistic sense.

Many technical processes (error dynamic, electromagnetism or structuralmechanic) can be simulated by commercially distributed software packets.These software packages can be coupled with the present method tooptimize in a probabilistic sense predetermined goals like aerodynamicefficiency, electromagnetic emission behavior or mechanical stability.

For many complicated processes no models exist. In these cases withRODEO a data based optimizing can be performed.

In the following possible applications for the present method are shown.These are merely examples. The actual application range for the presentmethod is much greater.

1. An airline wants to reduce its delays. Firstly possible influencefactors are determined like for example desk time, baggage dispatch,start slots etc. On many subsequently following days data of theseinfluence factors and of the resulting delays are collected. The presentmethod locates the greatest influence factors and performs a data basedoptimizing (see section “design of experiments”).2. The weight of a product should be minimized, the mechanical stabilityshould be not lower than a given limit. The wall thicknesses of theproduct do fluctuate, since the rolling machines merely guarantee acertain accuracy. Therewith the weight of a product also varies and themechanical stability merely can be provided with a predeterminedprobability. The mechanical stability can be calculated with afinite-element-package. The present method calculates the minimalexpectation value respectively minimal mean value of the weight (seesection “probabilistic design goals” item 2.3. Many technical apparatuses and measuring devices must fulfillpredetermined accuracy predeterminations. Many influence factors andtheir variability lead to the end accuracy. First the present methodlocates the most important influence factors in view of the endaccuracy. Second the variability (inaccuracies) of the single influencefactors can be maximized with the target, that the end accuracy keepsthe demanded value (see section “probabilistic design goals” item 6).4. The operating point of a plant should be determined. On the one handthe operating point should be optimal relating to one criterion, on theother hand the plant should be insensitive to fluctuations of theinfluence factors. The present method calculates the Pareto set out ofprobabilistic optimality and robustness of the plant. Basing on this theapplicant can decide by himself which compromise of optimality androbustness he elects (see section “probabilistic design goals” item 4).5. The crash behavior of a car is investigated. It is demanded that thenegative acceleration of a dummy does not exceed a certain value. Anessential influence factor is the sheet metal thickness, which is arandom variable because of the inaccuracy of the rolling machines. Nowit is a demand, that the expectancy value of the sheet metal thicknessis as small as possible, but the negative acceleration should not exceeda certain value with a pre-given probability. There exists aconventional method for the simulation of the crash behavior, theconventional method can be coupled with the present method (see section“probabilistic design goals” item 2).

A further example for a nonlinear technical system may be an antennaconfiguration, whereby an input parameter is the length of thetransmitter part and a target value is the transmitting power.

The present method is not limited by the application examples statedabove. The examples are merely seen as possible embodiments of thepresent inventive method.

The present method uses mathematical formulas, which are practicallyutilized, to improve all kind of nonlinear values of the nonlineartechnical system can be length, area, volume, angle, time, frequency,velocity, acceleration, mass, density, force, moment, work, energy,power, pressure, tension, viscosity, and all further physical kinds ofquantitiestechnical systems. Input parameters and/or target (see“Taschenbuch der Physik”, Kuchling, Verlag Harri Deutsch, Thun undFrankfurt/Main, 1985, chapter 3.6.). Examples for technical systems aretransport means like cars or airplanes, electronic circuits, powerstations, phones, turbines, antennas, fabrication processes of allindustrial goods and so on. In all cases input parameters and targetvalues are identified and used for optimizing. An improvement takesplace especially in comparison with conventional design.

An embodiment for the present method is a certain software named “RODEO”standing for “robust design optimizer”.

SUMMARY

According to the present invention for optimization of technical systemswith uncertainties an optimization model is proposed which uses a targetfunction comprising the expectancy value E(y) of the technical systemy=f({right arrow over (x)}) or the variance Var(y) or a combination ofboth values. A possible constraint can be used in the proposedoptimization model and can be a failure probability P_(f) being holdwithin a given tolerance. Expectancy value E(y) and Var(y) are given byformulas

E(y)=

f({right arrow over (x)})ρ({right arrow over (x)})d{right arrow over(x)}  (3)

Var(y)=

(f({right arrow over (x)})−E(y))² p({right arrow over (x)})d{right arrowover (x)}  (4).

The failure probability is given by formula

P _(f) =P(g({right arrow over (x)})≦0)=∫_(g({right arrow over (x)})≦0)p({right arrow over (x)})d{right arrow over (x)}  (6).

The integrals in equations (3), (4) and (6) usually can not beanalytically calculated. The overall method to solve the optimizationmodel can be merely efficient, in case the integral calculating methodsare efficient. Therefore used methods are subsequently described.

Methods for Calculating Expectancy Value and Variance

The used methods belong to the class of so called Response Surfacemethods. Specifically two variants can be used:

-   -   Taylor approximation    -   Hermite approximation

Using the taylor approximation the function f({right arrow over (x)}) issquarely developed:

$\begin{matrix}{{{f\left( \overset{\rightarrow}{x} \right)} \cong {f_{T}\left( \overset{\rightarrow}{x} \right)}} = {{f\left( \overset{\rightarrow}{\mu} \right)} + {{\nabla{f\left( \overset{\rightarrow}{\mu} \right)}^{T}}\left( {\overset{\rightarrow}{x} - \overset{\rightarrow}{\mu}} \right)} + {\frac{1}{2}\left( {\overset{\rightarrow}{x} - \overset{\rightarrow}{\mu}} \right)^{T}{\nabla^{2}{f\left( \overset{\rightarrow}{\mu} \right)}}\left( {\overset{\rightarrow}{x} - \overset{\rightarrow}{\mu}} \right)}}} & (7)\end{matrix}$

and the approximation f_(T)({right arrow over (x)}) is inserted intoequation (3) respectively (4):

E(y)≅

f _(T)({right arrow over (x)})ρ({right arrow over (x)})d{right arrowover (x)}  (8)

Var(y)≅∫(f _(T)({right arrow over (x)})−E(y))² p({right arrow over(x)})d{right arrow over (x)}  (9).

The integral in (8) and (9) can be exactly calculated now.

A higher approximation accuracy is achieved by a Hermite approximation.y=f({right arrow over (x)}) is approximated by a Hermite approach ofsecond order:

$\begin{matrix}{{f_{H}\left( \overset{\rightarrow}{x} \right)} = {{a_{0}{H_{0}\left( \overset{\rightarrow}{x} \right)}} + {\sum\limits_{i = 1}^{n}{a_{1i}{H_{1i}\left( \overset{\rightarrow}{x} \right)}}} + {\sum\limits_{i = 1}^{n}{\sum\limits_{j = 1}^{n}{a_{2{ij}}{{H_{2{ij}}\left( \overset{\rightarrow}{x} \right)}.}}}}}} & (10)\end{matrix}$

Herewith the Hermite polynomials are given by

H ₀({right arrow over (x)})=1,H _(1i)({right arrow over (x)})=x _(i) ,H_(2ij)({right arrow over (x)})=x _(i) x _(j)−δ_(ij).

By inserting of f_(H)({right arrow over (x)}) instead of f({right arrowover (x)}) in (3) respectively (4) an approximation of expectancy valuerespectively variance is achieved likely to the Taylor approximation.The coefficients in (10) are calculated for example by solving thefollowing least square problem:

$\begin{matrix}{\min\limits_{a_{0},a_{1i},a_{2{ij}}}{\sum\limits_{k = 1}^{p}\left( {{f_{H}\left( {\overset{\rightarrow}{x}}^{k} \right)} - {f\left( {\overset{\rightarrow}{x}}^{k} \right)}} \right)}} & (11)\end{matrix}$

for given interpolation points (evaluation places) {right arrow over(x)}^(k). Since functional evaluations are expensive an adaptive methodis used. Therewith firstly starting with only few interpolation pointsand adaptively adding further interpolation points as long as theapproximation of expectancy value respectively variance results inamendments. The adaptive optimizing method generates further advantagessuch that by the first optimizing steps the integrals can be calculatedwith a low accuracy but should be more accurate close to the solutionpoint. For high dimension problems the accuracy of the integralapproximation can be also adaptively fitted.

Method for Calculating the Failure Probability

For an easy representation it is assumed that the random variables areindependent and standard normal distributed. The methods are alsouseable with a common case. Merely a transformation must be executedbefore.

By a first step a point {right arrow over (x)}* of the highest failureprobability, a so called beta point is determined. This point isresulting from the solution of the following optimizing problem:

min ∥{right arrow over (x)}∥²

relating to

g({right arrow over (x)})≦0.  (12)

Two variants for calculating the integrals in (6) are used:

-   -   Linear approximation    -   Hermite approximation.

Being {right arrow over (x)}* the solution of (10) (beta point) andβ=∥{right arrow over (x)}∥. By the linear approximation (FORM: “firstorder reliability method”) g({right arrow over (x)}) is approximated by

g({right arrow over (x)})≅a ^(T)({right arrow over (x)}−{right arrowover (x)}).

Accordingly the following approximation for P_(f) is generated:

P _(f)≅(−β)

whereby φ is the distribution function of the standard normaldistribution.

A higher accuracy is achieved by using the Hermite approximationg_(H)({right arrow over (x)}) of the function g({right arrow over (x)})close to the beta point. Likely to the approach (10) the followingequation is achieved:

$\begin{matrix}{{g_{H}\left( \overset{\rightarrow}{x} \right)} = {{a_{0}{H_{0}\left( \overset{\rightarrow}{x} \right)}} + {\sum\limits_{i = 1}^{n}{a_{1i}{H_{1i}\left( \overset{\rightarrow}{x} \right)}}} + {\sum\limits_{i = 1}^{n}{\sum\limits_{j = 1}^{n}{a_{2{ij}}{{H_{2{ij}}\left( \overset{\rightarrow}{x} \right)}.}}}}}} & (13)\end{matrix}$

The coefficients are also determined by the least square problem (11).To evaluate the quality of the approximation the beta point {right arrowover (x)}*_(H) relating to the Hermite approximation is determined. Thebeta point {right arrow over (x)}*_(H) results from the solution of thefollowing optimizing problem:

min∥{right arrow over (x)}∥²

relating to

g _(H)({right arrow over (x)})≦0  (14).

Again an adaptive method is applied. Interpolation points are added aslong as {right arrow over (x)}*_(H) and the main curvatures in {rightarrow over (x)}*_(H) do amend. For evaluating the failure probabilityP_(f) the integral in (6) is transformed into:

$\begin{matrix}{P_{f} = {{\int_{{g{(\overset{\rightarrow}{x})}} \leq 0}{{\rho \left( \overset{\rightarrow}{x} \right)}{\overset{\rightarrow}{x}}}} = {\int{\Gamma_{g}\frac{\rho \left( \overset{\rightarrow}{x} \right)}{\rho \left( {\overset{\rightarrow}{x} - \overset{\rightarrow}{\beta}} \right)}{\rho \left( {\overset{\rightarrow}{x} - \overset{\rightarrow}{\beta}} \right)}{\overset{\rightarrow}{x}}}}}} & (15)\end{matrix}$

thereby Γ_(g) is the indicator function of g({right arrow over (x)}):

$\Gamma_{g} = \left\{ \begin{matrix}0 & {{g\left( \overset{\rightarrow}{x} \right)} > 0} \\1 & {{g\left( \overset{\rightarrow}{x} \right)} \leq 0}\end{matrix} \right.$

The approximation of P_(f) is

$\begin{matrix}{{P_{f} \cong {\int_{{g{(\overset{\rightarrow}{x})}} \leq 0}{{\rho \left( \overset{\rightarrow}{x} \right)}{\overset{\rightarrow}{x}}}}} = {\int{\Gamma_{g_{H}}\frac{\rho \left( \overset{\rightarrow}{x} \right)}{\rho \left( {\overset{\rightarrow}{x} - \overset{\rightarrow}{\beta}} \right)}{\rho \left( {\overset{\rightarrow}{x} - \overset{\rightarrow}{\beta}} \right)}{{\overset{\rightarrow}{x}}.}}}} & (16)\end{matrix}$

The integral in (16) can be efficiently calculated for example by“importance sampling”.

Monte Carlo methods are usable merely for small systems. For optimizingproblems they are not usable at all. “Standard” Response Surface methods(it means without adapting) are also still to large-scale, at least foroptimizing problems.

By using adaptive Response Surface Methods efficient methods also forcomplex optimization tasks with technical pertinent target functions(see above) and “chance constraints” (failure probability) beingconstraints are achieved. The adapting can be performed by two steps:

-   -   Adapting with the functional evaluation. The function is        adaptively interpolated (evaluated) at further points (places),        until a pre given accuracy of the target function respectively        of the constraint is achieved.    -   Adapting with the optimizing. At the first optimization steps        merely a low accuracy is necessary.

Using the adapting a given tolerance is achieved by minimal efforts.

1-13. (canceled)
 14. A method for optimizing a technical system that isa nonlinear technical product or process, comprising: describing thetechnical system by stochastic system output target values dependent onstochastic system input parameter values, thereby stating discretetechnical system dependencies; using Response Surface Methods based onthe discrete technical system dependencies to generate at least onecontinuous auxiliary function ({tilde over (f)}_(i)({right arrow over(x)}); {tilde over (g)}_(i)({right arrow over (x)})) for realdependencies of target values on the input parameter values (y=f({rightarrow over (x)}); z=g({right arrow over (x)})); using an auxiliaryfunction ({tilde over (f)}_(i)({right arrow over (x)}); {tilde over(g)}_(i)({right arrow over (x)})) to generate at least one optimizingparameter (E({tilde over (f)}_(i)({right arrow over (x)})); Var({tildeover (f)}_(i)({right arrow over (x)})); P({tilde over (g)}_(i)({rightarrow over (x)}))) optimized by an objective function, therebygenerating an additional discrete technical system dependence ({tildeover (f)}_(i)({right arrow over (x)}_(i)); {tilde over (g)}_(i)({rightarrow over (x)}_(i)); repeating both of said using operations toadaptively optimize the technical system, respectively additionallyusing the additional discrete technical system dependence ({tilde over(f)}_(i)({right arrow over (x)}_(i)); {tilde over (g)}_(i)({right arrowover (x)}_(i))), until a difference of successive optimized optimizingparameters (|E({tilde over (f)}_(i)({right arrow over(x)}_(i))−E(f_(i+1)({right arrow over (x)}_(i+1))|; |Var({tilde over(f)}_(i)({right arrow over (x)}_(i))−Var({tilde over (f)}_(i+1)({rightarrow over (x)}_(i+1))|; |P({tilde over (g)}_(i)({right arrow over(x)}_(i)))−P({tilde over (g)}_(i+1)({right arrow over (x)}_(i+1)))|) isas low as desired and a final additional discrete technical systemdependence ({tilde over (f)}_(p)({right arrow over (x)}_(p)); {tildeover (g)}_(p)({right arrow over (x)}_(p))) is useable as an optimaltechnical system operating point.
 15. A method according to claim 14,wherein the optimizing parameter includes at least one of a mean value(E({tilde over (f)}({right arrow over (x)}))) determined by E({tildeover (f)}({right arrow over (x)}))=

({tilde over (f)}({right arrow over (x)})ρ({right arrow over (x)}), anda variance (Var({tilde over (f)}({right arrow over (x)}))) determined byVar(y)=

({tilde over (f)}({right arrow over (x)})−E({tilde over (f)}({rightarrow over (x)}))²ρ({right arrow over (x)})d{right arrow over (x)},where ρ({right arrow over (x)}) is a probabilistic density function ofan input parameter distribution.
 16. A method according to claim 15,wherein the objective function is provided by minimizing or maximizingthe mean value (E({tilde over (f)}({right arrow over (x)}))) and/orminimizing of the variance (Var({tilde over (f)}({right arrow over(x)}))).
 17. A method according to claim 16, wherein the objectivefunction is alternatively or cumulatively provided by maximizingoptimizing parameters which are input parameter tolerances$\left( {\max\limits_{{{\sigma_{i,}j} = 1},\ldots,n}{\sum\limits_{i = l}^{n}\sigma_{i}}} \right)$and taking into account at least one constraint (P({tilde over(g)}({right arrow over (x)})≦0)).
 18. A method according to claim 16,wherein the objective function is alternatively or cumulatively providedby keeping an optimizing parameter which is a failure probability(P({tilde over (g)}({right arrow over (x)})≦0)) of at least one ofanother value, and the target value (P({tilde over (f)}({right arrowover (x)})≦0)), under or equal to a prescribed probability tolerance,with P({tilde over (g)}({right arrow over(x)})≦0)=∫_({tilde over (g)}({right arrow over (x)})≦0)ρ({right arrowover (x)}) and with ρ({right arrow over (x)}) a probabilistic densityfunction of the input parameter distribution.
 19. A method according toclaim 16, wherein the objective function is alternatively orcumulatively provided by minimizing a failure probability (P({tilde over(g)}({right arrow over (x)})≦0)) of at least one of another value andthe target value (P({tilde over (f)}({right arrow over (x)})≦0)) withP({tilde over (g)}({right arrow over(x)})≦0)=∫_({tilde over (g)}({right arrow over (x)})≦0) ρ({right arrowover (x)})d{right arrow over (x)} and with ρ({right arrow over (x)}) aprobabilistic density function of the input parameter distribution. 20.A method according to claim 19, wherein for mixed stochastic inputparameters ({right arrow over (x)}) and deterministic input parameters({right arrow over (x)}_(D)) a constraint depends on the deterministicinput parameters (h({right arrow over (x)}_(D))≦0).
 21. A methodaccording to claim 20, wherein the objective is provided by usingstochastic sensitivities of the optimizing parameters (E({tilde over(f)}({right arrow over (x)})); Var({tilde over (f)}({right arrow over(x)})); P({tilde over (g)}({right arrow over (x)}))).
 22. A methodaccording to claim 21, wherein normal distributed input parameters({right arrow over (x)}) form the probabilistic density function:${\rho \left( \overset{\rightarrow}{x} \right)} = {\frac{\sqrt{{Det}\Sigma}}{\left( \sqrt{2\pi} \right)^{n}}{\exp \left( {{- \frac{1}{2}}\left( {\overset{\rightarrow}{x} - \overset{\rightarrow}{\mu}} \right)^{T}{\Sigma^{- 1}\left( {\overset{\rightarrow}{x} - \overset{\rightarrow}{\mu}} \right)}} \right)}}$and/or beta distributed input parameters ({right arrow over (x)})comprise the probabilistic density function:${\rho \left( \overset{\rightarrow}{x} \right)} = {\prod\limits_{i = 1}^{n}{\frac{\Gamma \left( {\alpha_{i} + \beta_{i}} \right)}{{\Gamma \left( \alpha_{i} \right)}{\Gamma\beta}_{i}}{{x_{i}^{{ai} - 1}\left( {1 - x_{i}} \right)}^{\beta_{i - 1}}.}}}$23. A method according to claim 22, wherein the objective function isprovided by additionally or cumulatively using the cumulative densityfunction (

({tilde over (f)}({right arrow over (x)}d{right arrow over (x)})) of thetechnical system response.
 24. A method according to claim 23, furthercomprising generating the auxiliary function ({tilde over (f)}({rightarrow over (x)}); {tilde over (g)}({right arrow over (x)})) on based ondiscrete input parameter values ({right arrow over (x)}_(i)) anddiscrete target values (y_(i)=f({right arrow over (x)}_(i))).
 25. Amethod according to claim 24, wherein concerning changes of the inputparameters depending on time (Δ{right arrow over (x)}(t)) the inputparameters ({right arrow over (x)}) satisfy common stochasticaldifferential equations, whose density development is described by theFokker-Plank equation.
 26. A method according to claim 25, wherein theResponse Surface Methods are at least on of a Hermite, Laguerre,Legendre and Jacobi approximation.